Critical thermodynamics of the two-dimensional ±J Ising spin glass 
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We compute the exact partition function of 2d Ising spin glasses with binary couplings. In these 
systems, the ground state is highly degenerate and is separated from the first excited state by a gap of 
size 4J. Nevertheless, we find that the low temperature specific heat density scales as exp(— 2J/T), 
corresponding to an "effective" gap of size 2J; in addition, an associated cross-over length scale 
grows as exp(J/T). We justify these scalings via the degeneracy of the low lying excitations and by 
the way low energy domain walls proliferate in this model. 

PACS numbers: 75.10.Nr, 75.40.-s, 75.40.Mg 



Spin glasses 0,0 are strongly frustrated materials that 
have challenged statistical physicists since many years. 
In particular, there is still no consensus on the nature of 
these materials' phase diagram, a very basic issue. Sur- 
prisingly, open questions remain even in the case of two- 
dimensional spin glasses. For instance, there is a long- 
standing dispute [j, 0, concerning the ± J Ising spin 
glass: it is not clear what kind of singularity arises in its 
free energy at the critical temperature. 

In this work we reconsider the nature of these singular- 
ities using recently developed methods 0, 0] for comput- 
ing the exact partition function of square lattices with pe- 
riodic boundary conditions, focusing on the low T scaling 
properties of the model with binary couplings. We show 
that although the energy "quantum" of excitation above 
the ground state is 4J, such excitations behave as com- 
posite particles; in fact the specific heat near the critical 
point scales as if the elementary excitations were of en- 
ergy 2 J. We justify this picture using properties of exci- 
tations and domain walls in this model. Finally, the joint 
temperature and size dependence shows the presence of a 
characteristic temperature-dependent length that grows 
as exp(J/T), in agreement with hyperscaling. 

The model and our measurements - The Hamilto- 
nian of our two-dimensional (2c?) spin glass is 



Hj({o-i}) = -J2 J v a * a i 



(i) 



where the sum runs over all nearest neighbor pairs of 
Ising spins {pi = ±1) on a square lattice of volume V = 
L x L with periodic boundary conditions. The quenched 
random couplings Jij take the value ± J with probability 
1/2. The partition function at inverse temperature /3 = 
T- 1 is Zj = V r „ , e -l 3Hj ^' Ji ^ and can be written as 



Zj((3) 



Pj(X 



(2) 



Here Pj(X) is the polynomial whose coefficient of X p 
is the number of spin configurations of energy E = 



(-2L 2 + 2p)J. Saul and Kardar H |j| showed that de- 
termining Pj can be reduced to computing determinants 
which they did using exact arithmetic of arbitrarily large 
integers. More recently a more powerful approach has 
been developed 0, Q , based on the use of modular arith- 
metic to compute pfaffians. With this algorithm, one first 
finds the coefficients modulo a prime number, thereby 
avoiding costly arbitrary precision arithmetic. Then the 
computation is repeated for enough different primes to 
allow the reconstruction of the actual (huge) integer co- 
efficients using the Chinese remainder theorem. 

The algorithm proposed and implemented in 0, Q is 
powerful enough to solve samples with L « 100; the total 
CPU time needed to compute Zj grows approximately 
as L . In our study we have determined Zj for a large 
number of disorder samples at different lattice sizes: for 
instance we have 400000 samples at L = 6, 100000 at 
L = 10, 10000 at L = 30, 1000 at L = 40 and 300 at 
L = 50. The total computation time used is equivalent 
to about 40 years of a 1.2 GHz Pentium processor. For 
each sample we derive from Zj various thermodynamic 
quantities such as the free energy Fj(j3) = — /3 _1 lnZj, 
the internal energy Uj(/3) = (Hj), and the specific heat 
dUj/dT. We also study in detail the number of ground 
states and of excited states. Note that flipping any spin 
changes the energy by 0, ±4J or ±8J; the gap between 
the ground state and the first excited state is thus 4 J. 

Low temperature behavior of cy - The study of 2d 
Ising spin glasses has a long history. We will only discuss 
here results about the ± J model. It is generally agreed 
that this model is paramagnetic for T > 0, spin-glass or- 
dering arising only as T — > 0. The critical region thus 
corresponds to T — > T c = 0. Since there is an energy 
gap 4 J, the free energy should have a singularity of the 
form exp(— 4J/T). This is difficult to check, in particu- 
lar via Monte Carlo where the free energy is not directly 
measurable. Instead, it is better to concentrate on the 
specific heat density cy. For that observable, the differ- 
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ence between the models with bimodal (Jy = ±J) and 
continuous couplings is striking: in the first case cy goes 
to zero rapidly as T — > while in the second there is a 
clear linear behavior. 

Even though our computations provide us with the free 
energy, we also prefer to work with cy. Note that cy is 
related to a second derivative of the free energy so the 
corresponding singularities are directly related. Also, cy 
should provide a cleaner signal as irrelevant "constants" 
such as the ground state energy that fluctuate from sam- 
ple to sample have been subtracted out. Consider now 
any given sample. As T — > 0, we have the scaling 

( [H - (H)} 2 ) 16J 2 e 5l ~ s " e~ 4J / T 

where Sq and Si are the logarithms of the degeneracies 
of the ground state and first excited state energy levels 
for the given sample. (Sq and S± are microcanonical en- 
tropies; we have dropped the index J denoting a sample 
dependence.) Note that 4 J appears because it is the en- 
ergy gap in our system. It thus seems unavoidable that 
cy will have an exp(— 4J/T) singularity. Surprisingly, in 
1988, Wang and Swendsen [3j postulated that instead 

cy « T- p exp (-AJ/T) (4) 

with A — 2. They performed a Monte Carlo study in 
which A w 3 for most of the temperatures they could ac- 
cess, but their effective A drifted and their final predic- 
tion was A = 2 from an analogy with a one dimensional 
model (we shall come back to this later). This issue was 
taken up a few years later by Saul and Kardar 0, 0| who 
claimed A = 4; their work is based on exact computa- 
tions of partition functions and thus does not suffer at 
low T from the thermalization problems of the Monte 
Carlo approach. We are aware of no specific heat mea- 
surements in this model since. How could A not be 4? 
The subtlety is that we must take L — > oo at fixed T, and 
only after can we take T -> 0; indeed Eq. gj) assumes 
L = oo whereas Eq. JJJ) assumes T — > at fixed L. 

Using the algorithm in 0, Q , together with the avail- 
ability of cheap and powerful computers, we have ex- 
tended significantly the study of Saul and Kardar. For 
the sake of comparison, they had 80 samples at L = 20, 
22, and 24, and 4 samples at L — 32 and L = 36. (They 
also had samples for L < 18.) We go much beyond that, 
both in lattice sizes and in the number of samples we 
consider. In the left part of Fig. ^ we show our first anal- 
ysis of cy as follows. Let's set J = 1. When T — ► 0, if 
naive scaling {A = 4) holds, In (T 2 c v ) + 4/T w const, 
while In (T 2 cy) + 4/T - ^4 if A ^ ^ na i' ve and 
p = p n ai ve = 2. (The cy resulting from our exact parti- 
tion function computations has been averaged over dis- 
order samples.) In the plot we see that for any given 
lattice size, when T becomes small enough there is a sat- 
uration toward the naive scaling behavior, i.e., the points 
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FIG. 1: On the left: In (T 2 c v ) + 4/T versus l/T. On the 
right: -T\d(T 2 c v ) versus T. 

go to a constant value. The physically relevant regime is 
the thermodynamic limit, given by the envelope of these 
curves; this envelope does appear and seems to be linear 
in l/T. Note that the envelope emerges only on quite 
large lattices [L > 30); because of this, the true scaling 
escaped detection by Saul and Kardar. The straight line 
in the left part of Fig.^is our best linear fit to the L = 50 
data when (3 € [2.5,5.5]. It is a very satisfying fit and 
gives A — 2.02 ± 0.03, close to the integer value A = 2. 

We can also present the data in a slightly different 
fashion. In the right part of Fig.Q]we plot —Tin (T 2 cy) 
versus T. Here the coefficient A is given by the intercept 
of the envelope's extrapolation to T = 0, the left axis 
of the picture. We can distinguish three regions. The 
first region is for very low T values. Here the naive (non- 
thermodynamic scaling) with A — ^. na 'f V e — 4 is very 
clear. This region, where the intercept at T = is 4, 
shrinks to zero with increasing lattice size. In a second 
region we have the physical scaling; for the large lattice 
sizes we have, the value A w 2 emerges. Notice that this 
is the same region where in the left part of Fig. 2] the 
L = 50 data lie on a straight line. The third and last 
region corresponds to "high" T (T > 0.4) where one is 
far from the critical point and no scaling is apparent. 

Our conclusion here is that thanks to the larger sizes 
available to us and to a technique that does not suffer 
from low temperature critical slowing down, the thermo- 
dynamic scaling of cy is now finally clarified. 

Ground state properties - Our computations also 
give the ground state energies and degeneracies. Theo- 
retical arguments suggest that the mean ground state 
energy density eo has power corrections in l/T: 

e (L) = e* + a L- 2+&M . (5) 

We have = —1.4017(3) which agrees well with previ- 
ous work. We also find ©( e ) = -0.08(7); note that the 
prediction in is that 0^ e ^ = Odw, the exponent as- 
sociated with domain wall energies. Following the work 
of Hartmann and Young 9J, there is general agreement 
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that Odw — in the 2c? ± J model. Thus our estimate for 
is in excellent agreement with the conjecture in 

We have performed a similar study for the mean 
ground state entropy density so(L). We find Sq = 
0.0714(2) which compares well with the recent work 
of [13 in which s<$ = 0.0709(4). The fit also gives 
G)( s ) = 0.42(2), though if we take into account systematic 
effects we cannot rule out = 1/2. We believe that 
this large value, unrelated to 8dw, denotes the presence 
of a subtle organization of the ground states. 

Anomalous density of excitations The micro- 
canonical entropy S(E) of an energy level E is defined as 
the logarithm of the number of spin configurations having 
exactly that energy. Clearly, S(E) is obtained from the 
knowledge of Pj as computed in Eq. J3J). Of major inter- 
est is Si — So, the increase of entropy when going from the 
ground state energy E to the lowest excitation energy. 
In the pure ferromagnetic model, the lowest excitation 
corresponds to taking the ground state (all spins paral- 
lel) and flipping a single spin. This gives Si — So = In (V). 
One says that the excitations are "elementary" , and the 
system at low temperature is accurately described as a 
gas of independent excitations. 

The situation changes dramatically when a large 
enough fraction of J+j are negative, taking the system 
from a ferromagnetic to a spin-glass phase. In such a 
phase, the large V law of Si — So is modified. This is 
illustrated in Fig. [3 where we plot our numerical esti- 
mate of Si — So as a function of In (V). The dotted line 
is ln(y) while the dash-dotted one is 21n(V). We see 
that the true scaling behavior emerges only for large lat- 
tices and that the large V behavior is compatible with 
a 2 In (V) growth. How can one interpret this anoma- 
lous growth? Imagine classifying all excitations in terms 
of the size of the cluster of spins flipped when compar- 
ing to a given ground state configuration. (Naturally, 
one may ignore all clusters of zero excitation energy, and 
it is enough to focus on connected clusters.) Just as 
in the ferromagnetic case, some of the excitations corre- 
spond to single spin flips; there are 0(V) such objects. 
For any bounded-size cluster, the number of objects is 
0(V), leading to Si — So « ln(V). Since one has in- 
stead Si — So ~ 2 In (V) , finite-size clusters are irrele- 
vant: necessarily large scale excitations dominate the set 
of excitations of lowest energy. 

It is important to understand the nature of these 
large scale excitations, but unfortunately our computa- 
tional approach does not give us configurations, it merely 
counts their number. There are other ways to gain in- 
sight into this problem. To begin, we consider as in Q 
an analogy with the Id pure Ising model. In that system, 
when using periodic boundary conditions, the lowest ex- 
citation is composite, corresponding to a pair of kinks 
with a total energy 4J; however the "true" elementary 
excitations are single kinks, necessarily absent when us- 
ing periodic boundary conditions. It is easy to see that 
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FIG. 2: S(E + 4J) - S(E ) versus In (V) and the functions 
ln(V) and 2\n(V). 



for this Id model the quantity Si — So grows as 2 In (V), 
i.e., as in our 2d system. 

How may objects of energy 2 J appear in our 2d lattices 
with periodic boundary conditions for which the gap is 
4J? To answer this question, consider in a ground state 
configuration any connected cluster of spins and asso- 
ciate to its surface the corresponding closed path V on 
the dual lattice [llj. (The cluster's surface is the set of 
edges connecting the cluster to its complement.) When 
flipping the cluster, the change in the configuration's en- 
ergy comes only from those bonds crossing V\ in fact, 
for each such bond that is satisfied (JijSiSj = 1), the 
energy increases by 2 J, and otherwise it is decreased by 
2 J. It is easy to see that all clusters lead to V with 
an even number of bonds, and thus excitation energies 
are quantized in units of 4 J. However there are closed 
paths that are not associated with clusters: an exam- 
ple is a path that winds around one of the directions 
of the lattice! Such topologically non-trivial paths are 
called domain walls; when comparing periodic and anti- 
periodic boundary conditions, the set of bonds in the 
ground state that are changed from satisfied to unsatis- 
fied or vice-versa form exactly such a path. When L is 
odd, domain walls have energies ±2 J, ±6 J, . . . and the 
the quantum 2J appears. Of course, to have a physical 
excitation, one needs to introduce domain walls in pairs] 
then the flipped cluster of spins is topologically a strip 
with a surface in two pieces, one for each domain wall, 
while its energy is a multiple of 4J. Note that this is 
exactly what happened in the one dimensional case, the 
domain walls there being simply kinks which also arise 
in pairs. 

To justify the anomalous scaling of Si — So, we ap- 
peal to the fact realized a few years ago by Hartmann 
and Young |fj that low energy domain walls proliferate 
in the ±J spin glass. Let 5S be the typical entropy of 
a single domain wall; if we focus on those excitations of 
energy 4J associated with two domain walls of energy 
2J, we have an excess entropy AS = 2<5S. Our data 
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thus suggest 5S = In (L 2 ) ; this law can be interpreted by 
saying that SS is the sum of a In (L) term coming from 
the L possible mean transverse positions of the domain 
wall and of an additional In (L) term coming from the 
degeneracy (proliferation) at a given position. To extend 
this reasoning to the case of L even, we first remark that 
the domain walls there have energies 0, 4 J, ... To have a 
"strip" excitation, we need one domain wall of energy 4 J 
and one of energy. Undoubtedly, the entropy of these 
domain walls increases with their energy; a simple pat- 
tern is obtained if we conjecture that the excess entropy 
increases by In (L) every time the energy increases by 
the quantum 2 J. If this is so, the first domain wall will 
contribute 3 In (L) to the excess entropy and the second 
ln(L), leading again to the desired 21n(V) result. Such 
a conjecture is quite elegant and should be amenable to 
testing using a recent Monte Carlo method 01 • 

Finite-size scaling — Given the result for Si — So, we 
go back to Eq. (0) to understand the finite-size scaling of 
cy . When T — > and L — > oo simultaneously, standard 
arguments lead to 



T 2 c v (L,T) e 



2/3J 



T[L/k(T)\ 



(6) 



Here A(T) is a temperature dependent length that de- 
termines the cross-over between the thermodynamic scal- 
ing of cy (going as exp [—2/3 J] when L = oo) and the 
"naive" scaling as in Eq. ©. T is a finite-size scaling 
function; when its argument is large, L 3> A(T), we re- 
cover the thermodynamic limit and thus necessarily T 
must tend toward a constant. (Since cy is intensive, the 
L dependence must drop out.) On the contrary, when 
L <C A(T), we recover the behavior of Eq. J3J) where cy 
goes as exp (—4/3 J) but diverges as L — > oo. Interest- 
ingly, in this regime cy is not self-averaging and so one 
should apply finite-size scaling for the whole probabil- 
ity distribution of cy. Just as before where we computed 
S\ — So and not log(exp (Si — So)):, we focus on the typ- 
ical behavior and so we consider the median rather than 
the average of cy (this distinction is relevant only in the 




very low T, unphysical region, while it is irrelevant in our 
scaling region for T, say, close to 0.3). When using this 
data we have a very reasonable data collapse, consistent 
with Eq. © as shown in Fig. El We find that this median 
scales as L 2 exp(— 2/3J) at low T and thus !F(x) ~ x 2 as 
x — > 0. This then gives 



A(T) ~ exp(/3J) . 



(7) 



Summary and discussion — We have investigated the 
critical thermodynamics of the 2d Ising spin glass with 
binary couplings. Our main conclusion is that the specific 
heat density scales as cy ~ exp(— 2 J/3). This scaling is 
"anomalous" in the sense that it does not follow from the 
size of the energy gap (which is 4J). To find this scaling 
law, it is necessary to go to rather large systems, L > 30. 
We also found that the typical degeneracy of the first ex- 
cited level grows about L 4 times faster than that of the 
ground state level. We believe this high degeneracy has 
its roots in the proliferation of domain walls, two domain 
walls enabling one to define a composite excitation. Such 
a picture justifies the analogy with kink pairs proposed 
many years ago by Wang and Swendsen each do- 
main wall may indeed play the role of a kink, albeit with 
an additional entropy contribution. Finally, using finite- 
size scaling, we found a cross-over length scale A(T) that 
grows as exp( J/T); this divergence is exactly as expected 
from hyperscaling arguments. 

This work was supported in part by the European 
Community's Human Potential Programme under con- 
tracts HPRN-CT-2002-00307 (DYGLAGEMEM) and 
HPRN-CT-2002-00319 (STIPCO). 



FIG. 3: Data collapse plot of the finite size scaling function 
T [L/A(T)] with A(T) = exp{J/T) and J = 1. 
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